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ABSTRACT 


We discuss the direction dependent transfer of X-rays in a plane parallel 
atmosphere with a strong magnetic field perpendicular to the surface. We present 
a transfer formalism incorporating the full angular and polarization dependence 
of the cross sections, including vacuum polarization, for frequencies not too 
close to the cyclotron resonance. We treat the problem of a slab illuminated 
from below, and of a semi-infinite medium at constant temperature and density, 
and present numerical results for parameters typical of the hot polar caps of 
accreting magnetized neutron stars. Theoretical beam and X-ray pulse shapes 
are obtained for various models of X-ray pulsars, and the frequency and phase 
dependence of the pulse structure is briefly compared with observations 



I. INTRODUCTION 


The directionality of the radiation emitted from a magnetized atmosphere 

is an as yet not well understood problem, which has lately become of great interest 

because of its bearings on the beaming pattern of X-ray pulsars. With the advent 

of direct measurements of the magnetic field strength via cyclotron line features 

in several quasisteady, pulsating X-ray sources (Trllmper et al., 1978, Wheaton 

et al., 1980, Dennis et al., 1981), values of several x 10^^ Gauss have been now 

confirmed, which should strongly affect the transport properties of the atmospheres 

of these objects. Stimulated by these observations, some detailed radiative cross 

section calculations have also become available (Mesz^ros and Ventura, !978, 1979; 

Gnedin, Pavlov and Shibanov, 19/8; Ventura, 1979, Burner and Mfeszciros, 1979; 

Bussard, 1980; and Kirk and Mfesz^ros, 1980). These cross sections are strongly frequency 

and angle dependent. However, radiative transfer calculations used for models 

of X-ray pulsars have up until now concentrated on the frequency distribution of 

the radiation, using various angle averaging schemes (Ventura et al., 1978; Nagel, 

1980; Bonazzola et al., 1979; and Mesz^ros et al., 1980). The angular structure 

of a pencil beam was studied with approximate cross sections by Basko and Sunyaev 

(1975), through Monte Carlo simulations by Yahel (1980), and in the optically 

thick case by Kanno (1980). In the present paper we set up the transport equations 

for plane parallel finite and semi-infinite media, with full consideration of the 

directionality of the radiation. We have generalized Sobolev's (1963) integral 

equation treatment, to take into account the two normal polarizations of the 

medium and the angular anisotropy of the scattering and absorption cross sections, 

which is required for a realistic treatment of the magnetized neutron star 

atmosphere. Nagel (1981) has also attacked this problem by a different method, 

based on Feautrier's equations,, We have solved here the problem of a finite 

coherently scattering slab illuminated from behind, with internal sources of 

emission, and also the problem of the semi -infinite radiating medium, for a 
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magnetic field perpendicular to the surface (we discuss the case of the field 
parallel to the slab elsewhere). This approach is valid both for optically thick 
or thin situations and Ulows exact boundary conditions to be applied, namely, 
zero incoming diffuse radiation at the boundaries. In principle, by including 
thermal and quantum corrections in the scattering kernel, one could treat changes 
in frequency as well as direction in each scattering event. A full treatment is, 
however, bound to be numerically very difficult, and we have made here the simpli- 
fication of neglecting the frequency changes in order to concentrate on the inter- 
play of the polarization and direction changes. We have used the vacuum corrected 
cold plasma normal modes and assumed the electron distribution to be a one-dimensional 
Maxwellian. For photons near the core of the cyclotron resonance, e. g., Wasserman 
and Salpeter (1980), Kirk and Mesz^ros (1980), one would have to include both 
thermal and incoherence effects and also take into account departures from a 
Maxwellian distribution, e. g., Langer et al (1980). Our present calculation, 
therefore, does not apply for uj - However, over most of the continuum part of 

the spectrum, on which we here concentrate, we expect that the present simplifying 
assumptions should be fairly good. These calculations with the magnetic field 
perpendicular to the free surface provide a model for pencil beam patterns of 
X-ray pulsars in the continuum region of the spectrum. These beam patterns are then 
convolved with the rotation of the neutron star, assuming various aspects of the 
orientation of the magnetic and rotation axes to the line of sight, to give 
theoretical pulse profiles. An interesting result is that the contribution of the 
vacuum polarization to the anisotropy is quite strong, especially within ± 40 
percent of the cyclotron frequency. Also of interest is the appearance of double 
pulse structure, already indicated by Basko and Sunyaev (1975) on the basis of 
approximate calculations. Triple and double pulses can in fact appear, as the 
frequency becomes progressively less than the cyclotron frequency, and we discuss 
a possible method of determining the field strength B using this phenomenon. These 
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calculations also show the presence of phase dependent energy spectral distributions 
resembling observed ones. A large variety of pulse profiles can be generated with 
these very simple atmospheres by varying the density, the optical depth and the 
external illumination, so that detailed model fits may be ambiguous, especially 
if in the future one extends these calculations to inhomogeneous atmospheres. On 
the other hand, except for the observed asymnetry of some profiles (c.f.. Section V), 
it appears that these theoretical pulses provide a very flexible took for fitting 
many of the major features of observed pulses, 
n. NORMAL MODE RADIATIVE TRANSFER 

In the case of an isotropic, nonpolarized medium, the problem of the 
directionality of the radiation in a scattering atmosphere has been discussed in 
differential equation form by Chandrasekhar (1950), and in integral equation ^orm 
by Sobolev (1963). We extend here the integral equation fonnalism to a medium 
with two normal modes of propagation, with corresponding anisotropic and 

frequency dependent scattering and absorption cross sections, and a preferred 
direction given by the magnetic field. For the conditions in X-ray pulsars, 

throughout most of the continuin one has Uj+Uk)/2, where n^. , are 

the refractive indices and absorption coefficients respectively, e.g. Nagel (1980), 
MeszaVos et al. (1980). As shown bt Gnedin and Pavlov (1974), this reduces the 
radiation density matrix to diagonal form, so that one obtains a set of equations 
involving the normal intensities only. We shall assume a plane parallel con- 
figuration, with a magnetic field B directed perpendicular to the surface, along 
the z coordinate. In this case, the two normal intensities I, a>nd Ig (per unit 
circular frequency uj) depend only on z and G, the angle between the direction 
of observation and the field direction (z-axis). We can write 

cos fl cv, (s2)I, (z.ij) + ef(z,u) + 0 j*(z,q) + ej°(s2) 

dz 

cos P - ^2 (s2) Ig ( z,R) + e|(z,p) + e!,(z,S2) + ( q ) , 

dz • 


( 1 ) 
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where <v^(n) Is the extinction (scattering plus absorption) coefficient In cm"i 
for radiation of polarization 1 travelling In the solid angle dQ around n, 
and Is the emissivlty In erg/cnr“ sec'i ster for mode 1 Into that 

solid angle. There is an emissivlty term due to Internal sources (thermal 
emission), e°(q) » another due to external sources (outside Illumination) 
e*(z,f 2 ), and another due to scattering of diffuse radiation In the medium, 
e^(z,s2). We assume the thermal emission and absorption to be independent 
of z, which is valid for a homogeneous isothermal matter distribution. We find the 
Intensity scattered into d», c^(z,q) as follows. The energy of radiation 
of polarization j coming from ds? which Is stopped per unit volume and time 


Is (Yj(o')Ij(z,ii')dn’’. The energy scattered from mode j into mode 1 
Xji(o')cVj(n')Ij(z, n')dn', where 


is 
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We have then 


= f(7.fl) " »,')((!' )»i (ii' )I, (7. n' ) ■") ♦ 

w 


+ Xgi(Q' )l2(z, il' ) 


Xai(s/, n) 


( 4 ) 


= I, (7, 0' ) + 12(7, il ) " > 

dsi dQ 


For the emissivity due to external i lluinlnation, we denote the external 
radiation in the direction incident upon the free surface at z = 0 by 
I^(52 q). At a depth z, this gives rise to a diffused emission in the 
direction fi equivalent to 


e^(z. n) 


(“o> 


«) 


■o-i (52 q)z/C 0S Op 

dn 



( 6 ) 


Defining the source functions 


bJ (u) = ej (n)/fv.(Q), (6) 

Where the superscript n stands for s, * or o, we can rewrite the system 
(1) as 

cos fl dl^ (z t n)= -I. (z, il) + B. (z * ii) 
os (il) *'d“z 

cos 0 dl 2 (z, q)= -l 2 (z, il) ^ Bg(z, si). 


( 7 ) 



where 


B^(z, n) * b|(z, fi) B*(z, n) B°(n) 


( 8 ) 
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The formal solution of equation (7) is then 


1,(2. II) ■ IV, (n) ? 


COS 0’ 


(9) 


where z„, a constant determined from the boundary conditions. For these, 

we assume a slab configuration, of width z^, satisfying, 

0 


1^(0, sO = 0 for 0 < ri/g 

(10) 

(2p, ft) = 0 for 0 > n/g , 

that is, no incoming diffuse radiation at z * 0 and z = z ^ (!.| denotes the 
diffuse component only, the direct illumination from outside being I.). It 
follows that in equation (9), we must set 0 for 0 < n/2, and 2 ^ «z ^ 
for 0 > n/2. Replacing (9) into (4), (6) and (8) we obtain the integral 
equation for B^, 


B^(z, ft) 


= Pdft' 


V(z, ft) 

,rv. (ft) dft 


I 2 (z, ft' ) d(7p^ (ft| ftj 

~Cft) “^dsf 


+ B*(z, ft) + B°(ft) (11) 

Written out in full, with account of the boundary condition (10), this is 


B.(z,n ) = S / dep' 

’ 0 

r z 

. I J sin e'dfi 1 da..(i'i,ii) J (^, ((*i)d z 


f sin 6 d6 1 

V2 ■' dfJ 


do^.^iri,!^ 


7 ^ ^ COSG 


+ B* (z.ii) + B? (z,fti) 


( 12 ) 
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The external illuntlnation term, for a distribution of angles is, from (5), 
2n n/2 r j / . 

. .r ,r s.n . 

0 0 LI 

* "o , (13) 


while the thermal emission term is 


KTtTiT " rrTT;w/7^7CTr 


04) 


with P the Planck function for one polarization being 


P = (Zil)‘’c‘\nVCexp ( WkT)-l], 
w 


(15) 


where we have assumed implicity that the electron distribution is a one 
dimensional Maxwellian. The total intensity coming out at the other end 
of the slab (z = z^) is then, from equation (9) and the boundary conditions (10) 


“ (16) 


and the intensity reflected at z = 0 is 


AN ' n 


(17) 


cos e 


In (16) the first term is l^.(Zp,fi), the diffuse forward scattered intensity at 
Zq, and to make up the observed intensity, the second term representing the 
diminished incoming fTux has to be added. 
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ill. METHOD OF CALCULATION; FINITK SL AB 

A 

We Shan ffssume the Incoming flux i^(%) “ 
syinnetHc about the field B (z-direction). The integrated cross sections 
(e.g.. Appendix A) are fo-independant, crj (n) * while the differential 

cross sections da^ j(n' ri)/dn depend on From (13), B*(z,q) will therefore 

not depend on and from the iterative construction of B^(z,q) implicit in (12), 
neither will it, B^(z,a) * B^(z,o). The <•'' integral in (12) acts therefore only 
on dr^ j(fj'n)/dfi, leading to 


2 1 

B^(z,u) = £ / du' 

’ j=l 0 


^ 6,n^,u) /° } 


+ B* (z,u) + B° (u) 


(18) 


where u = cos o, we have defined 


G,j(u' .11) 


2n 

r df(T 


1 dcrjj(u^ ,u ) 1 




LU 

dfi 




<jpij(ulu ) 
du' "" 


( 19 ) 


and for the 6-integrals between n/2 and n we have performed the change of 
variable u' -♦ -u^ and used the synmetry properties (e.g.. Appendix A) 
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G^j(l/ ,u) * ,u) 

Pf^(-u) » cT^(u) (20) 


o(^(-u) « 'V^(u) 


The corresponding B 


* 

1 


Is 


1 

B*(z,u) « dUp[Gi^(M^. 

+ Gb^(Uo, 




( 21 ) 


We have solved numerically the system of equations (18), under the boundary 
condition of a given incoming beam I^(Uq) impinging on the slab surface at 
z = 0 and obtained the reflected (z = 0,li/2< 0 < n) and outcoming (z=z^,0<o ,<ii/2) 
beams lj(u), lj(u), by setting a grid of z, u and u' values, and solving iteratively. 
Some care must be given to the treatment of the exponential, c.f. Appendix B. The 
nature of the cross sections, as determined in particular by the value of 
the ratio of the frequency to the cyclotron frequency, influences the accuracy 
of the calculation for any finite grid. For very large optical depths, especially 
if the cross sections are very different for the two modes, as is the case for 
o)/u^ « 1, the finite slab numerical solution becomes progressively less accurate, 
and a different set of equations is more advantageous, as described in the next 
section. 
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IV. SEMI-INFINITE MEDIUM 


Let us first establish the form of the equilibrium solution at very 
great optical depths, wj z » 1, o^z» 1. At these depths, the Influence 
of B* on the solution should be negligible. We notice that the exponential 
terms will be non-negllgible only for z-z' uV<v^(u^). If the optical depth 
Is large for all u, this means we may set B(/, u' ) ^ B(z, u' ) and take It 
out from under the Integral sign in (18). This 1s further justified by the 
consideration that physically one expects B not to depend on z at great 
depths. If we take z at some arbitrary point, we may set the upper and 
lower limits of / as + ® and write 


B. M = } du' fe, (U',) G., ,u) ; I"' d/ t 

+ %(.') eE,(u',u) Je-|z-2'l"a(u')/u' caip-dz'] 

««QO 

Bs(u) • } du' ["b, (u-.)Gis(u' ,u) r aiilll dz' ♦ 

0 ^ -00 U 

+ B,(u',) %u(u'.u) ; =6^^ dz-]. 


( 22 ) 


that is, 

Bi(u) = 2 } Bi(uf ) Gij(u^,u) du" + 2j Bg (u" )(^i (uf ,u)du" 

0 0 

1 1 

Bg(u) = 2 / Bi (u" ){^j (u". ,u) du" + 2/ BgCu" 

0 0 


( 23 ) 



n 


We seek the conditions for an Isotropic solution, B^(u) =1^* constant. 
Because of the sytmnetry properties, (u.u'ydu' * dcr^^ (u',u)/du and 
dcTy (u,uO/du' = d«Tj^(u\u)/du (i/j),(c.f. Appendix A). We can write 


1 dtr (i,^,uM TT- i ^^,»{u,uM ./ 

” B| 2 ” " "y "1 I. 1 i ■ ^ ■. 1 du EU 2* " ■ X& m mm . .i‘<. 

^ ‘ 0 du' 0 '^TTiiT du» 

du + 2.1 


(24) 


^ ° ^ 5Tu7 


From the definition 


a.(u) = /du' . 

^ du' .[ du' 


i 


and the fact that cb. .(u,u' )/du' = cb . . (u,-u' )/du' , so that J du' 

1 J I J Q 

these integrals, we obtain 

^ r ] d(T, ,(u,u' ) ^ ^ ) ^ / 

•^i f ^ ^2 f —^^dlT ° 

‘ 1 da, Ju.u' ) ‘ I da (‘J.'J') 


(25) 


= ^Jdu' in 




<^TuT 


du' - 0 


(26) 


In equilibrium therefore we have 


Bj - , Ij = la 


( 27 ) 


as one might expect. 


Let us seek now approximate expressions for B^(u,z), valid at large 
optical depths, but not quite as large that it would lead to (27). We 
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consider a semi -infinite medium, extending from z = -«> to z * +Zp. The 
equations are 

1 Zo 

B,(2,u) =J'((u' G,.(u'.u) ; dz' + 

' n ' nn 1 1' 


(28) 


♦ Jdii' Gj,(b',u) / dz' 2z^B^(z',u')e-l"-^l"'2<'''>/“' 


At z = -«» Bj = Bj =62=62, and we may choose z^ large enough that at z= 0 
this is still true. This will be the case if for all u, we have 


Zq » Max(i/o’^(u)), i = 1,2 . 


(29) 


Normalizing to = 1, (28) becomes 


1 r ^0 

B.(z,u) = Jdu' Gi .(i/,u) f Bj(z' 

n r\ 


.u')e-l^-2'I“i(‘'')/u' i^dz' + 


0 

z 


+ /du'^Gi .(u',u)e-"i^“'5^/''' + ^.(u\u)e-'^<"')^/" 


yl 


(30) 


The last tei'm, which comes from the integral between z = -» and z = z^ under 
replacement of z^ by +®, plays the role of boundary condition. 


V. RESULTS 


We have performed calculations for a magnetic field strength of B = 0.1 B 

c 

= 4„4 X 101® Gauss, perpendicular to the surface of the atmosphere, for 
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finite slabs and for semi -infinite media, with grids in 2 , u and u' typically 
of 30 points each. The densities used reflect the values expected In an 
accreting X-ray pulsar, varying from 10*^^ to 10** cm-3. However, we are not 
aiming here at providing a model for any specific source, but rather have in 
view a study of the properties of the radiative transfer in a general situation 
involving a magnetic field. 

An interesting question is the effect of a magnetized atmosphere in 
reprocessing the directional distribution of an incident beam of light of 
arbitrary angular distribution. In Figure 1 we have plotted the shape of 
the outgoing flux at a photon energy of 10 keV as a function of angle for 

aa , 

a slab of density n = 10 cm" , T = 10 keV and Thomson optical depth tj = 7. 

The curve labelled p is the outgoing intensity for an incoming beam distri- 
bution which is a delta function in angle (actually a step function of width 15°), 
centered at n = 0° (parallel to the field). The curve labelled "cos" is for an 
input varying as cos n, while curve A is for an isotropic input. The curve B is 
the outgoing intensity for zero input (i.e., due only the self-emission of the slab, 
breinsstrahlung modified by scattering). The curve C is for a semi-infinite 
medium of the same density and temperature. The general effect, of course, 
is towards a broadening of the incoming beam, which increases with the optical 

depth of the atmosphere. For an input delta function at 0 = 0° an optical 

_aB a 

depth of = 7 (based on = 6.6 X lO” cm , not on 't(g)) is however not yet 
enough to smear out the incoming beam entirely. Similarly, a cosine input 
is modified but higher optical depths would be required to reproduce an 
isotropic input, A semi -infinite medium produces a beam similar in shape 
to the isotropic input in a slab. A self emitting slab (D) produces a hollow 
beam, as already pointed out by Basko and Sunyaev (197b) » which reflects the 
structure of the cross sections. If one convolves these beam shapes with the 
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rotation of the neutron star, assuming that both poles emit at the same 
relative intensity, we obtain the pulse shapes of Figure 2. We have used 
here as illustration two different choices of aspect angles ij (between the 
line of sight and the spin axis) and ig (between the spin axis and t), namely 
(il » I 2 ) =(60°, 46°) and (75°, 45°). The effective angle at which the 
beam is sampled by the observer is given by 

cos 0 = cos ij cos ig + sin ij sin ig cos '•n, (31) 

where is the azimuthal angle around B (O when the beam is pointing closest 
to the line of sight), and if n > n/2 one takes n -♦ ii-n. The profiles are in- 
variant under interchange of ij and ij. 

In Figures 3, 4 and 5 we show the outcoming fluxes in both polarizations 
at photon energies of 1, 10, 20, 30 and 70 keV. The cyclotron frequency 
is near 50 keV for our field value of B = 0.1 B^, and in the present 
paper we avoided this frequency range, since to study it properly would 

require if,troducing a number of complicating factors, e.g., Kirk and M^sz&ros 
(1980). Already for 30 and 70 keV, if the photon propagation angle, during 

its random walk, comes close to 0 or n , the photon would have entered the 
core as defined by Wasserman and Sal peter (1980), so that within the present 
approximations frequencies closer to than these could lead to errors. The 
density and temperature were taken as n = 10^-^ cm"^ and T = 10 keV. Figure 3 
shows the flux from a slab of Tj = 7 when illuminated by an isotropic black- 
body flux (T(^^ = 5 keV). Figure 4 is the flux from the same slab when the 
incoming luminosity is zero (self-radiation only). Figure 5 is the flux from 
a semi-infi ni te medium of the same characteristics. These fluxes were computed 
using the cross sections of Appendix A, the full lines being calculations 
including vacuum polarization, and the dashed ones neglecting it. One sees 



that for an external illumination exceeding the luminosity of self- 
radiation (Figures) the hollow beam structure giving rise to double pulses 
(visible in Figure 4) is not present, except to some degree at such frequencies 
where the background intensity (Tj^j^ = 5 keV) is less than the plasma emission 
(T = 10 keV), as is the case at 70 keV for polarization 1 (Figure 3a). It 
is interesting to compare our Figure 4 with Nagel's (1981) slab results with 
comparable parameters (his Figure 7a). Some difference is inherent in the 
fact that we included vacuum polarization, which in Nagel's calculation was 
set equal to zero, and also his profiles are summed over polarizations and 
normalized to the blackbody intensity at each frequency. Allowing for this, 
the two methods seem to lead to fairly similar results. It is also interesting 
to compare Figures 4 (slab of Ty = 7) and Figure 5 (semi-infinite medium), 
both without external illumination. As one would expect, the central hollows 
in the beams tend to fill up as the path length increases. 

It is useful also to present beam shapes normalized to unity, for 

singling out the purely angular features, and in Figure 6 we have done this 

for a semi-infinite medium of n = 8x10^3 cm"^, B = 0.2 Be, T = 10 keV for 

various values of , where in this case is 100 Kev. The full lines 

are calculations including vacuum polarization, and dashed lines are without. 

One sees that for a given density and temperature, as one decreases -o/uj the 

H 

single lobe evolves toward a triple lobe, and as one decreases even 
further it becomes eventually a double lobe. Kanno (1980) has computed lobe 
shapes for the same parameters, using an approximate analytic procedure, which 
may be compared directly with our Figure 6. It is interesting to see that 
his approach leads to qualitatively similar predictions concerning the pulse 
multiplicity* but there are quantitative discrepancies. For instance for 
values of ('n/dt^) of 20 and 30, with vacuum polarization, Kanno predicts 
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the peak of the side lobes at about 45° for bothi and normalized intensities 
of about 0.16 and 0.18 respectively. Our calculation on the other hand pre- 
dicts 38° and 32° with intensities of 0.6 and 0.78 respectively. In general, 
his method produces rather narrower central lobes than ours, e.g., 8° at 
intensity 0.5 vs. our 19° at 0.5 for = 10, and much flatter and less 

intense sidelobes. A qualitative difference is that we find normalized 
central lobes which widen with increasing , while Kanno's method 

gives the opposite behaviour, at least for = 10,20,30. 

While the nonnalized lobes reveal the angular changes with frequency, 
the intensity information contained in Figures 3, 4 and 5 is linked with 
these. The effect of vacuum polarization in the normalized lobe plot of 
Figure 6 is to narrow the angular distribution and lower the intensity. In 
the non-renormalized lobes of Figures 3, 4 and 5 we see however that the 
inclusion of vacuum polarization reshuffles the energy between the frequencies 
and the polarizations. This is because the radiative cross sections in- 
cluding emission are strongly affected, both in their angular and frequency 
behaviour by the inclusion of the magnetic vacuum polarization, e.g. Mesz3ros 
and Ventura 1978, 1979. In the usual plasmas, in the magnetoionic, cold 
plasma limit one has negative electrons gyrating all in the same direction 
and giving rise to predominantly elliptically polarized modes, so that only 
one mode is resonant, that for which t gyrates in the same sense as the e". 

In the magnetic vacuum of an X-ray pulsar, however, the virtual e^-e" pairs 
can become more important than the real e"o One has now charges of both 
signs, almost equal in number and gyrating in opposite directions. This has 
the effect of driving the normal modes towards becoming nearly linearly 
polarized, and it also makes the second (ordinary) mode resonant, so that the 
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relative importance of both modes is changed. The cross sections, as seen 
in appendix A, are direct functions of the normal polarization vectors 
and reflect therefore these changes. Although, as discussed for Figure 6, 
the general effect of the vacuum is to narrow the beam, this appears so be- 
cause in this plot the flux is normalized. However the inclusion of the 
vacuum also alters the relative intensity at a particular frequency, be- 
cause of its effect on the mean free paths and thermal ization lengths, etc. 
This can lead to an increase of intensity at some angles and a decrease at 
others, e.g. for 0<20° and G>20° at 10 keV in Figure 5a, or only to an in- 
crease, e.g. 20 keV in tite same graph, or only a decrease, e.g. at 70 keV. 
This naturally shows itself also in the pulse profiles that we discuss 
further below. 

We have plotted in Figure 7 the differential photon energy density Uyj 
inside the atmosphere as a function of depth z (linear scale), for the 
inch'vidual polarizations and for their sum, at a photon energy Ruj = 10 keV. 
Curves "a" corresoond to the slab with background of Figure 3, curves "c" 
to the slab without background of Figure 4, and curves "b" to the semi- 
infinite case of Figure 5. The depth z runs from zero to z^, and in the slab 
cases the radiation escapees from both ends, but in "a" there is an external 
illumination incident at z = 0. The slabs have = 7, while for the semi- 
infinite case "b", the value of z^ was chosen so as to satisfy the 
condition (29). For instance, at 10 keV, we took z^ = 1400 cm or = 94. 

The number of iterations used for solving the integral equations has to be 
chosen large enough to ensure that information is propagated across a 
distance equal to the slab length, or in the semi-infini te medium across a 
thermali zation length. Thus at 10 keV, where a 22 >>a]i'uji 2^021 >>K 2 >>ici, for 
the slab this implies > (ctj 2 2 ^)' . Typically we used 40 to 60 iterations 
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for the slab models, and over 400 for the semi -infinite models. 

We then calculated pulse profiles, shown in Figures 8, 9, and 10 corresponding 
to the intrinsic beam shapes of Figures 3, 4, and 5, by convolving them with the 
rotation of the neutron star at different values of the inclination angles ii andi;> 
of equation (31), indicated at the top of each figure. The total fluxes (1 + 2 .) 
at each frequency (70, 30, 20, 10 and 1 keV from top to bottom) were all normalized 
to unity, to emphasize the difference in shape. The full lines are the calculations 
including vacuum polarization, while the dashed lines are the corresponding cal- 
culations without vacuum polarization. It is seen that at most frequencies > 10 keV, 
the inclusion of vacuum polarization introduces significant differences, and the 
more so the closer one is to the cyclotron resonance (50 keV for these figures). 
Depending on the aspect angles ii, ig, and the frequency, the pulses can appear 

broader or narrower than for the zero vacuum case. For the high intensity cases 
(slab with background. Figure 8 , and semi-infinite case. Figure 10) at 10 keV the 
vacuum narrows the pulse, but in the low intensity, self-emitting case (Figure 9) 
it broadens it. It also changes the depth of the hollow, where present. At 20 
and 30 keV, the vacuum in general narrows the pulses, but at some angles this effect 
reverses. At frequencies larger than the cyclotronic, in general it seems to 
broaden them, but in practice this is more difficult to observe since the total 
photon number in the actual observations is lower, and the noise higher. These 
changes introduced by the vacuum polari zation* reflect the corresponding ones in the 
cross sections, discussed previously. The cross sections can become larger at 
some angles at the expense of other angles, but since the transfer involves a 
complicated interplay of photon creation, scattering and polarization exchange, 
the relative intensities in both polarizations at different frequencies are changed, 
so that depending on the viewing angles i^ and i 2 and the phase angle this can lead 
to either narrowing or broadening of the pulses. 

We have also calculated pulse profiles based on the intrinsic lobes of 
Figure 6, with n = 8 x lO^s cm"3, B = 0.2 Be (co^ = 100 keV), T = 10 keV, shown 
in Figure 11. The aspect angles ii, ig used were 40/33 and 78/86, and the energy 
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variable was taken to be u)/u)^|, for comparison with Kanno (1980). Our calculation 
predicts wider central peaks, and the transition to a double peak occurs at lower 
frequencies than In his approach. 

From the point of view of the application to specific models, one sees from an 
examination of the sample of theoretical pulses described above that a broad range 
of pulse shapes can be reproduced by changing a few basic parameters of even the 
simplest, uniform atmospheres. One feature observed in many X-ray pulsars is an 
asymmetry of the pulse shape around phase 0, or 1, (Pravdo et al., (1977), Pravdo 
et al., (1979), and this does not appear In our theoretical pulses,. One way to 
reproduce such asyiiwetries might be to have an inhomogeneous flow along the azimuth 
of the accretion column. Another possibility woulri be an asyi?rietr1c(non-dipole) 
magnetic field. Without, however, going into the asyiwnetry, which is not very 
pronounced in many pulsars, we can attempt a qualitative comparison with some 
known objects. Two well studied objects are A0535 + 26 (Bradt et al„, 1976) and 
3Ud900-4Q (McClintock et al., 1976). The first shows no Interpulse, the main pulse 
being single at higher energies and splitting into two at lower energies. Looking 
through Figures 8, 9, 10 and 11 we see that the self-emitting slab of Figure 9 
around 50/20, 30/30, or the semi-infinite medium of Figure 10 at 30/30 or of 
Figure 11a around 40/33 shows this behavior. In 3U0900-40, at higher energies 
there are single main and interpulses, of not too different height, splitting 
at lower energies into double structures, possibly going through a multiple pulse 
stage in between. This can be seen in the theoretical profiles of Figures 9, self- 
emitting slab at 80/60 or 75/45, and the semi-infini te medium of Figure 10 at 80/60 
or of Figure 11b at 78/86. Kanno in fact suggested, on the basis of his approximate 
calculations, semi-infinite atmospheres at 40/33 and 78/86 as possible candidates. 

It is important to notice that one has the choice of several models, involving 
different physics., A semi-infinite medium kept at 10 keV is of course an ideali- 
zation in an accreting neutron star. In principle, the heated region should have a 
depth comparable to the stopping length of incoming protons, which depends on a 
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number of other uncertain variables (whether they are stopped by particle collisions* 
or If there Is a shock, etc.) but perhaps 5-50 gm/o'n^ is the depth of the actual X-ray 
atmosphere. If the heating Is more or less uniform throughout, a self-emitting 
slab could then be used as a first approximation. On the other hand. If the 
heating occurs mostly towards the end of the stopping length, a slab with a black- 
body Illumination from below might serve as approximate model. If downward thermal 
conduction were very effective, which again Is uncertain, a semi-infinite atmosphere 
might be appropriate. Better (and needless to say more difficult) models would have 
to consider also density, temperature and field gradients, as well as departures 
from plane parallel symmetry. It Is encouraging however that even these very 
simplified Models are able already to account for several major features of the 
observed pulses. 

Of potentially great observational Interest are the frequency of onset of 
the multipulse structure, and the phase dependence of the energy spectrum. The 
first of these could, in principle, serve to determine the magnetic field strength, 
if one did not know it from, say, a cyclotron line measurement. The reason is that, 
for a given density, temperature and geometry, the multipulse behavior sets in at a 
particular value of u)/u)|^. This method is not applicable however to all pulsars, 
since the aspect angles have to be within a certain range (these are aspect angles 
at which no multipulse structure shows up)a It also makes a difference whether 
slab or semi-infinite models are used. However, if one is able to adjust the total 
luminosity, by making reasonable guesses for n , t and Z^, a fit of a particular 
sequence of pulses at different energies with a theoretical sequence at different 

would determine u)^, or B. The other effect is that the spectrum observed varies 
with phase, as a look at Figures 3, 4 and 5 shows. This effect is known to occur 
observatlonally, e.g., Pravdo et al., (1977), a hardening of the spectrum occurring 
towards the middle of the pulse of HER X-1. Without going Into details, which would 
require a number of further considerations, one sees that for instance the self- 
emitting slab of Figure 4 shows a definite hardening towards 9=0, and the semi- 
infinite model of Figure 5 also has a similar tendency. 
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To suninarize, we have Investigated here sonie of the main properties of the 
basic radiative transfer In a pencil beam situation. These results give one a 
feel for the elements one should use In a realistic model of an accretion column. 
Insofar as the structure of the radiation field. It Is worth noting that the 
amount of computation required In this approach Is significantly less than In 
Monte Carlo photon escape calculations. The results discussed here need to be 
extended to include the fan beam case (field paralled to the free surface, c.f. 
also Nagel 1981), and to take account of Incoherent scattering, a calculation now 
in preparation by the present authors. It is clear. In conclusion, that the 
peculiarities of the cross sections in a magnetic field, as well as the structure 
of the atmosphere, play a large role In determining the pulse shapes of accreting 
neutron stars. Detailed models of astrophysical sources, incorporating some of the 
elements described above, should be possible in the near future. 
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APPENDIX A; CROSS SECTIONS 


We list here the expressions for the cross sections necessary to arrive 
at the symmetry relationships used In the text (for a discussion, see Misziros 

and Ventura 1979, Ventura 1979). In this appendix, u is no longer cos 0 as in 
the text, standing for a different quantity, c.f. (A2). 

The rotating components of the extraordinary (1) and ordinary (2) 
polarization eigenvectors ^ , P are 

P = (ei^•^ eV*“, ei_.») 

e »2 “ -C fv sin = -C sin o 

(A1) 

e»^ = C 2‘^ e'-'to cos « + 1], e* = C 2"^ e'^'^ Ccos e - o'] 
el ='c 2"^' cos 0 - Ij, e* = C 2"^ e'*’^'''' [cos n + cv] 

where 0 is the angle between the magnetic field and the wave vector, '•o is 

the axial angle around B, and 

^ fv(0) = b(l-[l + C = (1 + ^2)4^ 

b “ b(0) = 2"\sin*0/cos 0)u^[l + 36(l-u)/(uv)], (A2) 


U = V = («)p/(o)®. 


where and are cyclotron and plasma frequencies. The 8- parameter 
is the vacuum polarization correction, of value 
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« . (45n)’‘ (e*/Fc)(B/B^)* ~ 0.5 X 10‘*(B/B^)*. (A3) 


where B^ is the critical field B^ » m“cVeli - 4.4 X lO^’G at which the cyclotron 
energy equals me*. 


The differential scattering cross section from mode i at angle o into 
mode j at angle f>' is 


dn' ° 




J' 


1 + u- 


1- U‘ 


(A4) 


♦ 2 a,j(b,j *c,j)cosP t (b,j ♦c,j)“cos“pt(b, .-c, .)Mn“Pl 


where p = -rn, r^ = e*/mc®, 



leil 

lef 1 




i4i 

lej'l 

(1 + uh"‘ 

(A5) 


ki 

l-i'l 

(1 - U^)’‘ . 



the prime indicating that e' must be taken. The '•''-integrated differential 
cross section is 


da. .(0,e' ) 
sin 0'" di^ 


nr^fZ a. .* . (b.^ . c.^)* . (b.^ - c^^)*} = 


e^ e^l 


i 


= Zrrr^CdeMle'J I)*•^( 


1 u- 


(A6) 
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The integrated polarization exchange cross section is 


..,(8) ■ J" stn d«- • 

0 sin 0' d0' 

^ " (Hu?)« (1-u^ 

where A;J, = (3/4) / d(cos 0 ' ) |e)J(0\i-'' ® 0)|", and 


(A/) 


“ 1 

r A‘* = 1 

" ^cv ' • 


j=l 


(AS) 


The integrated cross section sunined over outgoing polarizations is 


"J ° 'j, + = <T^Uei: 


J la 


J |3 






(A9) 


The free- free absorption cross sections, in the limit where one sets 
Gaunt factors equal to unity (an approximation which should not introduce 
much greater deviations than some of the other ones made), can be taken as 
(e.g., Ventura, Nagel and M^szJiros 1979). 


Ki(0) = a.(0)(3ir/2)(z e^/Tiv^^)(ngC3/!«3)(l-e‘^'’/'^^), 


(AlO) 


where v^^ = (nkT/2m)‘' and n^ is electron density. 
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APPENDIX B; SKETCH OF THE NUHERICAL_ME mOIL 

The integrals in (18) require some attention to the fact that the kernel 
K^(u) = (cv^ (u'‘)/u- )e“^i becomes unwieldy when <y -♦ «. Let us define 

B+(J.L) = KB(J + 1. LM B(J + 1,1)] (3,) 

B_(J.L) - KB(J + 1. L) - B(J,L)]. 


where B(J, L) is the value of B computed on the point J, L of the integration 
grid. We have terms of the type 


/du' G(u'. u) / Mz'. „-)e-(^-z'Mu')/u' 


NL-1 

^ &= 0 L'k 


(L' + l)k J-1 J'h 

r du' G(u' ,u) S r B(z* ,u' ) 

j"=i (y-i)h 


(B2) 



where we have treated for simplicity one polarization only, and dropped the 
index i. Here h = z^/lNJ-l) is the z step size, and k is the u = cos o 

step size, and J = z/h. We may assume B and G to be sufficiently smooth that 
to a good approximation the integral (Bg) may be written 


1 2 NL-1 J-1 

/ du'.../ dz'.. = T r E [[G(L'k.u) + G((L' + l)k, u)] 

0 0 L'=0 J'=l 


[B+(J',L + 1) + B^(J'.L)] 


(I/+ l)k J'h 

J du' J 

L'k (J'-l)h 



The double integral becomes 


(L'+l)k J'h (L'+l)k 

r du' r g-(z-z')cv(u')V^2., P 

V ^ . V y ^ t-i 

L'k (J'-l)h L'k 


[-g-(z-h(J'-l)_g-(z-hJ')j^ 
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and defining the functions 




kL 

r 

1c(L. 


1 ) 




0^ * 


(B4) 


the integral (B3) becomes 


1 z NL-i J-1 

rdu'.. f dz'.. = i S E f[G(L'k.u) + C.((L>l)k,u)] 

0 ‘ 0 L'=l J'=1 


[B^(J',L+1) + B^(J'.L)] (B5) 


There is another similar integral, with z' < z' < z^, which analogously 
becomes 


1 ^0 NL-1 NJ-1 

J du\. J dz\. - i- r i: {[G(kL'.u) + G(k(L'+l).u)] 

0 z L'=l J"-1 


[B+(J\L+D + B^(JM)3 [.ji+(J'-J+l,L) + (B5) 


The functions can be computed to the required accuracy and stored at 

the beginning of the computation. For slabs of Ty < 3 this method is accurate 
to within a few percent. 

When the condition z^o'(u) » 1 is encountered (too wide a slab, or semi- 
infinite medium), the above method is no longer adequate. To deal with this 
case, instead of refining the grid size (which in fact does not help much), it 
is more advantageous to take into account the effects of the first derivative 
of B. The integrals (B4), (B5) can be written then as 
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r du\ . r dz' „ . , 


NL-1 J-1 A 


T: {[B^(y ,lM + M (J' ,L' 

L'=i y=i 


U^iO-y ,L')->i-^(J-J^l,L')] 

A A 

+ B_(J\L' )[.!.+ ( J'-J+1,LM-^.(J'-J.L']]G(L\u), (B6) 

A A 

where B+(J,L) = +• 8^(Jt L + 1)], G(J,L) = ^[G^(J,L)+G(J,L+'I)], and 

♦ =■ f"'. , ™ e-^C-XJ-^Wu 


Similarly, 


'M(k-l) 


' 0 NL-1 NJ-1 A A 

/ du' .. J dz\. = 7; E f[B^(y ,L' ) + {y y )] 

0 z L'=1 J'=J 


[>j)^(j'-J + 1. lM - + 2, lMI- 


- B_(J^ ,L' )]]G(L' ,u) 


From the definition of tp_, B^ and B_ one sees that for cv(u)Z » 1 the 
integral (B6) becomes 

1 z NL-1 

;du'... Jdz'.. = k t: KB(0. l') + B(J. + 1)] G(L^ u), (B8) 

0 0 L^-1 

1 z„ NL-1 

^du^.. Tdz'... = k i; L') + B(J, L^ + 1)] G(L', u) (B9) 

0 ^ L'=l 

These integrals are the same as one would have in the asymptotic case of 
equation (23), i.e., 

1 

B(0. u) = /G(u' ,u)B(J, u')du' + B*(J, u), (BIO) 

0 

with an added boundary condition term B*. 
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FIGURE CAPTIONS 

Figure 1 Outcoming beam shapes as a function of angle 9 for atmospheres of 
n = 10^2 cm-3, T = 10 keV, C: semi-infinite atmosphere, B: self-emitting 

of Tj = 7 and zero external illumination. A: same with external illumination, 

isotropically distributed in angle, at blackbody temperature = 5 keV. 

"Cos": same if external illumination is distributed as cos 0. "6": same 

if it is distributed as 6(0), la) Polarization 1* extraordinary, lb) Polari- 
zation 2, ordinary. The flux is per unit circular frequency (d, that is erg 
cm"2 sec-i sec. The flux is per unit circular frequency w, that is erg cm-2 
sec"^ sec. 

Figure 2 Pulse shapes obtained convolving the beam shapes of Figure 1 with 
rotation of the neutron star, 2a) for inclination angles ii, ia = 60°, 45°, 

2b) for ii ia = 75°, 45°, 


Figure 3 Outcoming beam shapes as a function of angle 0 for a slab of = 7, 
n = 1Q23 cm"3, T = 10 keV, B = 0,1 B^ = 4.4 x 10^2 ggygs (hw^ 50 keV), photon 
energies 1, 10, 20, 30 and 70 keV, illuminated by an incoming isotropic beam of 
blackbody radiation at T^^^ = 5 keV. 3a) polarization 1. 3b) Polarization 2. 

The flux is in erg cm"2 sec-^ sec, the full lines are calculations including 
vacuum polarization, and the dashed ones are what one would obtain if this is 
neglected. 

Figure 4 Same as Figure 3, slab of Tj - 7, but with zero external illumination, 
self-emission only. 


Figure 5 Same as Figure 3 but for a semi -infinite medium. 
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FIGURE CAPTIONS (Cont'd.) 

Figure 6 Outcoming flux from a semi -infinite medium of n » 8 x cm"^, 

Ts.10 keV and B = 0.2 Be TOO keV), normalized to unity. Full lines 

include vacuum, dashed ones do not. The label on the curves is now (w/u)^)"*. 

Figure 7 Photon energy densities u^[erg/cm^s"^] as a function of depth z at 
frequency 10 keV, for the slab of Figure 3 (curve a), the slab of Figure 4 
(curve c) and the semi-infinite medium of Figure 5 (curve b). Polarizations 
1, 2 and their sum 1 + 2 are given. 

Figure 8 Pulse shapes (polarization 1 and 2 summed) normalized to unit flux 
for 70, 30, 20, 10 and 1 keV (top to bottom) obtained convolving the beam shape 
of Figure 3 (slab Tj = 7 with background illumination Tj^j^ = 5 keV) with the 
rotation of the star, at different values of ij, ia, labelled at the top. 

Full curves: including vacuum polarization. Dashed lines: without it. 

Figure 9 Same as Figure 7 for the self-emitting slab of Figure 4. 

Figure 10 Same as Figure 7 but for the semi -infinite medium of Figure 5, 

Figure 11 Pulse shapes normalized to unit flux for the semi -infinite medium of 
Figure 6, for (w/w^)"^ = 10, 20, 30, 40 and 100 from top to bottom, and aspect 
angles of 40/33 and 78/86. 
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